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We report on a study of the superfluid state of spin-polarized atomic ®Li conflned in a magnetic 
trap. Density profiles of this degenerate Fermi gas, and the spatial distribution of the BCS order 
parameter are calculated in the local density approximation. The critical temperature is determined 
as a function of the number of particles in the trap. Furthermore we consider the mechanical stability 
of an interacting two-component Fermi gas, both in the case of attractive and repulsive interatomic 
interactions. For spin-polarized ®Li we also calculate the decay rate of the gas, and show that within 
the mechanically stable regime of phase space, the lifetime is long enough to perform experiments 
on the gas below and above the critical temperature if a bias magnetic field of about 5 T is applied. 
Moreover, we propose that a measurement of the decay rate of the system might signal the presence 
of the superfluid state. 

PACS numbers: 03.75.Fi, 67.40.-w, 32.80.Pj, 42.50.Vk 



I. INTRODUCTION 

One of the most important objectives in the study of 
dilute gases, has been the achievement of Bose-Einstein 
condensation (BEG) in bosonic systems. Indeed, decades 
of experimental research finally lead two years ago to the 
observation of BEG in three different systems of alkali 
metal gases ^■^Rb, ■^Li and ^SNa This success has 

triggered a large amount of interest in the field of ultra- 
cold atomic gases. Although the study of properties of 
these degenerate atomic Bose gases is vigorously being 
pursued at the moment, trapping and cooling of Fermi 
gases might also provide new and exciting physics. In- 
deed, in a previous theoretical study we showed that a gas 
of spin-polarized atomic ^Li becomes superfluid at den- 
sities and temperatures comparable with those at which 
the Bose-Einstein experiments are performed . 

This superfluid phase transition, which is similar to the 
BGS transition in a superconductor, occurs at such high 
temperatures due to the fact that ^Li has an anomalously 
large and negative (triplet) s-wave scattering length a 
This scattering length is a measure for the interatomic 
interactions, and its sign implies that this interaction is 
effectively attractive, which is a first requirement for a 
BGS transition to occur. For other atomic species, the 
transition temperature is in general very low, because of 
the fact that the scattering length is of the order of the 
range of the interaction ry and the dilutcncss of the gas 
requires that the Fermi wavenumber kp <^ l/^'y- So 
for example in the case of deuterium, it was concluded 
already some time ago that the observation of a BGS 
transition is experimentally impossible 

The ^Li atom has nuclear spin i = 1, and electron spin 
s = 1/2. Gonsequently the atom has six hyperfine states 
|1) to 1 6), for which the level splitting in a magnetic field 
is shown in Fig. |l|. The upper three levels |4) to |6) can 
be trapped in a static magnetic trap, whereas the lowest 



three hyperfine levels prefer high magnetic fields and are 
expelled from a magnetic field minimum. 
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FIG. 1. Energy of the six hyperflne states of Li in units 
of the hyperfine constant a^/, as a function of the magnetic 
fleld. 

The simplest way to create a degenerate Fermi gas is 
to trap just one low-field seeking hyperfine state, and for 
the sake of stability of the gas, the doubly polarized state 
|6) — \ms — 1/2, TOj — 1) is most suitable. However, due 
to the Pauli exclusion principle, two fermions in the same 
hyperfine state can interact with each other at best via 
p-waves, and if this interaction is effectively attractive, 
the onset of the formation of Gooper pairs occurs at a 
temperature 



Tc~ — exp 



2(fcF|a|)3. 

where ep — h^kp/2m is the Fermi energy of the atomic 



1 



gas, and a the p-wave scattering length. For ^Li this p- 
wave scattering length of the triplet potential is approxi- 
mately — 35ao, where ao is the Bohr radius, and kpa ^ 1 
in general. As a result the critical temperature for such 
a doubly spin polarized ^Li gas is extremely low. At 
present, a reasonable number for the density of trapped 
atomic gases is lO^^cm"^, leading to ep/kB — 600nK, 
and fecial ~ 7 x 10"'^. The corresponding critical tem- 
perature is clearly unattainable. 

In the case that more than one state is trapped, Bara- 
nov et al. [Q| predicted a considerable increase in the 
above (p-wave) critical temperature as a result of the fact 
that two atoms in the same spin state can now also at- 
tract each other through the exchange of a phonon (den- 
sity fluctuation) in another hyperfine state. It was found 
that in this case the transition temperature 



T.c^l^exp -13f^ 



where a now corresponds to the s-wave scattering length 
for collisions between the two hyperfine states. Never- 
theless, using again a density of 10^^ cm~^ for each spin 
state and the value a = — 2160ao for ^Li we find that 
kF\a\ — 0.43, and it is easily verified that also in this case 
the critical temperature is out of reach experimentally. 

Therefore, the most promising approach is to consider 
a Cooper pair of two atoms in different hyperfine states, 
since then s-wave pairing is allowed. In this case |q| 



Tc ~ — exp 

Kb I zfci?|a| 



resulting in a much higher critical temperature than in 
the previous two cases. In particular, we envision to trap 
^Li atoms in the states |6) and |5). Experimentally, this 
might be achieved most easily by first trapping only one 
hyperfine level, and then applying a noisy rf-pulse to cre- 
ate an incoherent mixture of atoms occupying these two 
hyperfine levels Note that this situation has in fact 
already been realized in recent experiments with *^Rb 
atoms, although using a different technique ||l0|| . 

In a recent publication Modawi and Leggett propose 
to trap ^Li atoms in three instead of two hyperfine states 
pT| . The advantage in such a system is that the effect of 
fluctuations is reduced somewhat, but the disadvantage 
of trapping more hyperfine states, is that the number of 
channels by which the gas can decay increases consider- 
ably. There are not only more possibilities for two-body 
collision processes, in which one or two electron spins are 
flipped and the corresponding atoms are expelled from 
the trap, but also three-body recombination processes 
are now no longer suppressed. Therefore, at present, it 
seems to be most favorable to trap only two hyperfine 
states, and the most suitable candidates are the states 
|6) = |m, = 1/2; m, = 1), and |5) ~ \m, = 1/2; = 0), 
because for this combination the decay processes due to 
two-body interactions can be suppressed most. The ap- 
proximate sign in the last expression indicates that in 



the state |5) there is for ^ ahj a small admixture 
with the spin state \ms = —l/2;mi = 1). Although this 
admixture can be neglected for most purposes, we come 
back to its importance for the stability of the gas shortly. 

As explained above, in a two-component spin-polarized 
atomic ^Li gas. Cooper pairing will occur only between 
atoms in different spinstates, while there is almost no in- 
teraction between two atoms in the same spin states. For 
notational simplicity, we also refer to these states as | "f) 
and I I), and the densities of atoms in these two hyper- 
fine states are denoted by and , respectively. Notice 
that since the two states are electron spin-polarized, the 
strength of the interatomic interaction is indeed char- 
acterized by the s-wave scattering length of the triplet 
potential Vt{v), and it is exactly this number which is 
anomalously large and negative in the case of ^Li. 

The aim of the present publication is threefold. First, 
the homogeneous calculation of Rcf. Q needs some im- 
provement, due to the fact that the interatomic interac- 
tion potential for ^Li has recently been determined more 
accurately The most up-to-date value of the s-wave 
scattering length is a = — 2160ao, where oq is the Bohr 
radius. This change in a not only affects the critical 
temperature but also the decay rates of the gas. Second, 
we want to take the effect of the inhomogeneity of the 
trapped gas into account, and in the third place, we look 
for a signature that signals the presence of the superfluid 
phase in the gas. 

The paper is organized as follows 



II A 



we 



In Sec 

consider in some detail the decay processes limiting the 
lifetime of the gas. Subsequently, we briefly summarize 



the theory for the homogeneous Fermi gas in Sec. [IB 



and improve the results obtained earlier for the critical 
temperature, using the mos t up-to-date interatomic po- 
tential for ^Li. In Sec. II C we consider the mechanical 

In particu- 



stability of a weakly interacting Fermi gas. 
lar we also consider a gas with positive s-wave scattering 
length, and show that in the unstable part of the phase 
diagram, a spinodal decomposition can restore the sta- 
bility of the gas in this case. 

In future experiments the atoms are likely to be 
trapped in an external potential that roughly has 
the shape of an isotropic harmonic oscillator V{r) = 
imw^r^, and which causes the gas cloud to be inhomoge- 
neous. Therefore, the last part of this paper is devoted to 
the study of an inhomogeneous two-component Fermi gas 
at and below the critical temperature, and in particular 
we will again concentrate on ^Li. The numerical calcula- 
tions will be performed in the local density approxima- 
tion, which is valid if the correlation length f over which 
the particles influence each other is much smaller than 
the typical trap size I = y^Ti/mu! over which the density 
of the gas changes. A similar calculation for the nonin- 
teracting case has been performed recently by Butts and 
Rokhsar |l2|. In addition, the case of purely repulsive in- 
teractions has been studied by Oliva in the same way in 
the context of possible experiments with spin-polarized 
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atomic deuterium ||T^. In Sec. [II A, we briefly repeat 
the i ngred ients for the local density approximation. In 



Sec. [II B we calculate the critical temperature of the 



func tion of the number of trapped atoms, and 
in Sec. [II C we study the gas in the superfluid state. 
Density profiles for the gas as well as for the BCS order 
parameter are presented. In Sec. IV we devote a discus- 
sion to the issue of how to detect the superfluid phase 
and to distinguish it from the normal phase. We end the 
paper with a summary of the main conclusions. 



II. HOMOGENEOUS FERMI GAS 

We first consider a homogeneous, dilute gas of (elec- 
tron) spin-polarized ^Li atoms. Since the gas is dilute, 
the atoms will interact with each other mainly through 
two-body collisions. These two-body collisions can be 
represented on the mean-field level by a local potential 
with a strength given by the two-body scattering matrix 
rp2B _ /^Tiafi^ jjji^ where m is the mass of the particles 
and a is the scattering length The sign of a de- 

termines whether the two-body interaction is effectively 
repulsive (a > 0), or attractive (a < 0). 

Before going to a description of the gas in the nor- 
mal and superfluid state, we consider an aspect that is 
experimentally of some importance, namely the lifetime 
of the gas. The large s-wave scattering length has on 
the one hand the advantage of having many thermalizing 
collisions between the particles which is required for ef- 
ficient evaporative cooling, but on the other hand there 
will also be relatively many inelastic collisions which can 
cause spin- flips within the atoms. If the electron spin 
of an atom is inverted, the atom will be lost from the 
trap, and consequently these inelastic processes limit the 
lifetime of the gas. In the next subsection we explain in 
more detail which decay processes dominate in a mixture 
of ^Li atoms in the hyperfine states |6) and |5). 



A. Decay rates 

The electron spin and nuclear spin quantum numbers 
of the two trapped hyperfine levels for 3> auf are 
given by 

|6) = |m, - 1/2; m, = 1) 

|5) - |m, = 1/2; m, = 0) +0+|m, = -1/2; m, = 1), 

where 9^ ~ ahf / {2^/2^eB) is inversely proportional to 
the applied magnetic field _B, so for sufhciently strong 
magnetic fields the admixture of |5) with the high- field 
seeking part is small and the gas can considered to be 
electron spin-polarized. For such large magnetic fields, 
the energies of these two hyperfine levels are given by 
^& — OLhf l'^, + HeB and £5 ~ fieB, respectively. 



Since the two atoms in state |5) and |6) will interact at 
the low temperatures of interest solely via s-wave scatter- 
ing, implying that the spatial part of the two-body wave 
function is symmetric under the exchange of atoms, the 
spin part of the wave function must be anti-symmetric, 
i.e. 



|{6,5}_) 



5>-|5)|6)] 



1 

111; 11) + 61+100; 22), 



(1) 



where in the last line we used the basis \SMs; I Mi) with 
S = Si -f S2 and I = ii -I- 12 the total electron and nu- 
clear spin of the two colliding atoms, and Ms and M/ 
the corresponding magnetic quantum numbers along the 
direction of the magnetic field. 

The decay rates for the transition from the state 
l^m, {a,/3}) with orbital quantum numbers I and m to 
a state |/'m', {a', /?'}) with quantum numbers V and m', 
is essentially given by Fermi's Golden Rule and results in 
the expression ||l^ 

I 1 2 

\Tl'm'{a' ,l3'},lm{a,l3}{Pa' ,P' , (2) 

for the zero-temperature limit of the rate constant for this 
process. Here Ti,jn'{a' ,i3'},hn{a,i3}iPa' ,f3' ,0) is the two- 
body scattering matrix at zero energy such that the in- 
coming particles have zero relative momentum, and the 
magnitude of the relative momentum of the scattered 
particles is Pa',0>- 

As was explained for example in Ref. Q , there are ba- 
sically two ways in which collisions cause the atoms to be 
lost from the trap. First of all, the central (singlet and 
triplet) interaction V = Vs{r)V'-^'^ +FT(r)'p(^) induces 
transitions between different hyperfine levels. Since this 
interaction cannot change the total electron or nuclear 
spin angular momentum, and the hyperfine level |5) has 
a small admixture with the state |rn.s = —1/2, to,; — 1), 
only transitions |{6,5}_) |{6, 1}-), where |1) ~ jm^ = 
— l/2,TOi = 1) —6^\ms = 1/2, TOi = 0) are allowed. Sim- 
ilar to Eq. (|l|), the total spin state |{6, 1}-) is given by 



|{6,1}_) = |00; 22) -0+111; 11) 



(3) 



Combining Eqs. ([^) and (||), we find that the spin part of 
the transition matrix Ioo{6.i},oo{6.5}(P6i, 0) contributes a 
factor 9^ times the exchange potential V^^^ (r) — Vr (r) — 
Vs{r), i.e. the difference between the triplet and singlet 
potential. To calculate the spatial part, we must use for 
the relative in- and outcoming scattering wave functions 
with orbital quantum numbers I and to and total electron 
spin S the normalized expression 



Vis 



(r) 



i%m{r), 



(4) 



where ipl^"^ (r) denotes the in- and outcoming solutions to 
the radial Schrodinger equation with the singlet or triplet 
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interaction. Using furthermore that the relative momen- 
tum pqi after scattering is due to the energy difference 
^6 ^ ^1 — 2fieB, we find thatpgi = v^2m/vS. Combin- 
ing all expressions into Eq. (|2|), we obtain that the rate 
constant due to exchange interactions is given by 

'(-l(r,p6i)|FT(r)-Fs(r)|vI/W(r,0)) 



, ^ n 



3*2 



'■hf 



^U{r,pei)\V'^^r)\¥+l{r,0)) 



(5) 



In Fig. ||this exchange rate as a function of the magnetic 
field is shown (curve 1). 

The second way in which collisions cause decay of the 
gas is due to magnetic dipolar interactions V^. As will be 
shown, of the various dipolar interactions, the contribu- 
tion due to electron-electron dipolar interactions is most 
important. For this dipolar interaction, we have Mm 




AMs 



2-AM, 



■Ami: 



AMs' 



(6) 

where the tensor operator Sj^^^j^j^ can be thought of as 
arising from the coupling between Si/h and S2/h, the 
Pauli spin matrices describing the electron spin of the 
two atoms, to a tensor of rank 2. For the scattering state 
|{6,5}_) ~ |11;11), the dipolar interaction can change 
the (total) electron spin projection Ms of the two atoms 
by an amount AMs = ~1 for a one spin-flip (Isf), or 
AMs = —2 for a two spin-flip (2sf) process. Therefore, 
the one (two) spin-flip dipolar interaction couples the in- 
coming wave function with approximately S = 1, Ms = 1 
to the final state having S — 1, Ms — {Ms = — 1). As a 
result, the outgoing wave function is in the state |10; 11) 
for one spin-flip, and in the total spin-state |1 — 1; 11) af- 
ter the two spin-flip interaction. The Clebsch-Gordon 
coefficients for each process are given by -^/S/IO and 
•\/3/5, i.e. the spin part of the transition matrix con- 
tributing to the decay rate is a factor of \/2 larger for 
the two spin-flip process than for the one spin-flip pro- 
cess. Moreover, the energy released in a one spin-flip 
process is only half of the energy released in a two spin- 
flip process. Therefore we find that pl^/f^, — \/2mp,eB 

whereas p^^/p, = \JAm^cB. We thus arrive at the conve- 
nient relation that G'^'f {B) = 2G^'f {2B) and it suffices 
to calculate only the one spin-flip decay rate. 

Performing a similar calculation as in the case of the 
exchange decay rates, the one spin-flip decay rate be- 
comes 



G^'f {B) = 2TT^n'^my/2mfieB 



12 



*{,+)(r))(10;ll|Ef_i|ll;ll) 



— V2TO/ieB — (rao) , 



(7) 



where 



{r2of = 



dr 



is the radial electron-electron dipolar element. In Fig. 
the one spin-flip decay rate constant is shown as curve 2. 

At this point it can be understood that the electron- 
electron dipolar interaction gives the largest contribu- 
tion to the dipolar decay rates. Decay due to the 
electron-nucleon interaction occurs for example via the 
|{6, 5}_) |{6, 4}_)-channel. However, the correspond- 
ing decay rates are smaller by a factor of {^j^ / ^eY' — 
20 X 10~® and thus completely negligible. This also im- 
plies that a mixture of |6) and |5) atoms cannot achieve 
equilibrium in the spin degrees of freedom within the life- 
time of the gas. This is completely analogous to the re- 
cent experiments with two condensates in different spin- 
states performed by the JILA group [pl)[. 




B(T) 

FIG. 2. Decay rate constants due to exchange (curve 1) and 
one spin-flip processes (curve 2) as a function of the applied 
magnetic bias field. 

Fig. ^ shows that the lifetime of the gas is of the order 



of 1 s for a density 



10 cm and a mag- 



netic bias field of 5 T. Although this would provide ample 
time to perform an experiment, a much shorter lifetime 
may be adequate. For successful experiments we not only 
have to require that the time between thermalizing col- 
lisions is small compared to the lifetime of the gas, but 
also that the time scale for formation of the Cooper pairs 
obeys this condition. The latter is anticipated to be of 
0{?i/kBTc) and therefore in our case much longer than 
the time between collisions. Nevertheless, for a density 
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JT-s = "-6 — 10^^ cm~^, we have that Tc ~ 11 nK, and 
h/kBTc is only about 0.7 ms. Our calculation overesti- 
mates the spin-exchange rate constant below a magnetic 
field of 1 T, and we estimate that a bias field of 0.2 T 
would give a lifetime of at least 1 ms. 

In the next subsection we consider the microscopic the- 
ory that describes the Fermi gas in the normal and the 
superfluid state. We only apply the BCS theory here. 
The influence of fluctuations j|] will be addressed in a 
future publication. 

B. BCS transition 

For the homogeneous case, and taking only two-body 
interactions between atoms in different hyperfine states 
into account, the gas is described by the following hamil- 
tonian 

= ^ / dx Vi(x) Ma ) V'c<(x) 



where tIq, = (■0ct(x)-0Q(x)) is the equilibrium value of 
the density of atoms in state |a), and equivalently Aq = 
Va{i! \ (x)-0t (x)) is the equilibrium value of the BCS order 
parameter |16|] . The chemical potential of each hyperfine 
state has now been renormalized to [i!^ = [i^ — T^^n^a 
to include, on the mean-field level, all two-body scatter- 
ing processes with particles in state | — a). The factor 
the two-body scattering matrix, and 
has been substituted for Vq to incorporate correctly all 
two body processes into the calculation. Note that the 
same substitution should not be performed in the expres- 
sion for Aq , because all two-body interactions are already 
going to be included by the BCS-treatment as we will see 
below ||l^ . Due to the nonequilibrium in the spin degrees 
of freedom, both chemical potentials /i'^ and fi'^ need not 
be equal, and therefore the densities of atoms in the re- 
spective hyperfine level can be varied independently. 
Substituting for the operator ipl^ the expression 

^t(x) = ^^a^e-*-, (10) 



J dx' l/T(x-x')Vi(x)V^„(x')V-a(x')^a(x). 



(8) 



In this expression, "f and J, refer again to the two hy- 
perfine states involved. The field operators '0q(x) and 
tpai^) obey the usual Fermi anti-commutation relations, 
and denote the annihilation and creation of a fermion at 
position x in hyperfine state \a) with chemical potential 
fia- The interparticle potential can be approximated by 
a local potential, Vt(x — x') ~ VoS{x — x'), where the 
constant Vq is a measure of the strength of the interac- 
tion. We come back to the precise value of Vb shortly, but 
it is in any case negative to account for the effectively at- 
tractive nature of the triplet interaction. The integration 
over x' in the hamiltonian is then trivial. The next step 
in a mean-field treatment of the hamiltonian in Eq. (|^), 
is to develop the operator products ipl,ipa and ■(/'aV'-a 
around their mean values by substituting 



and 



To first order in the fluctuations, we are left with the 
following effective mean-field hamiltonian 



H = 




- A^a ^a(x) 



+AoVj(x)^|(x) -f A*i/;i(x)i^|(x) 



Vo 



(9) 



where aj^ ^ creates one particle in spin state \a) with mo- 
mentum ?ik, the hamiltonian in Eq. (^) becomes 



H 



A* 



Ao 



ak,T 



— n^niT 



(11) 



where ek — h^\i^ /2m is the free particle energy of a par- 
ticle with momentum ?ik. The density of atoms in state 
I a) is determined by 



(12) 



Since the effective mean-field hamiltonian in terms of the 
operators aj^. ^ and ai^a is non-diagonal, one cannot di- 
rectly calculate the expectation value (aj(. c,ak,a)- 

This is, as usual, resolved by first applying a Bogoli- 
ubov transformation according to pq] 



k,i 



-Wk&k,T + ""k^-k 



(13a) 



(13b) 



to diagonalize the hamiltonian in Eq. (|l l|) . After per- 
forming this unitary transformation, we require that the 
hamiltonian in terms of the new quasiparticle operators 
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6k, T and 6^j. only has diagonal elements, and further- 
more that these operators again obey the usual anti- 
commutation relations for annihilation and creation op- 
erators. This determines the values of the yet unknown, 
and in principle complex, constants Wk and Wk- The latter 
constraint requires that the constants Uk and Uk must sat- 
isfy the relations | Uk P + K'k P = 1 , and the requirement of 
diagonality of the hamiltonian after the transformation 
leads to the condition |ukP = ^(1 + Ck/VCk + l^op), 
introducing ^k — ^k — e_Fj i-e- the free particle energy 
relative to the average Fermi level tp = (/x'^ -I- 

The eigenvalues corresponding to the Bogoliubov quasi 
particles are then given by 



|Ao| 



(14) 



where nia — ±1/2 for a =1? i: respectively. Furthermore 
Sep = fJ-'^ — fJ-'^ is the difference in Fermi levels of the two 
hyperfine states. The dispersion relations of Eq. ( p^ are 
depicted in Fig. ^ for equal (Figs. ^ and |^) and unequal 
densities (Figs. ^ and |^(i) both with zero (Figs. |^ and 
^) and nonzero Aq (Figs. ^6 and^), respectively fl^ ]. 
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FIG. 3. Bogoliubov dispersion hLJu,^a for a) ni = ni, 
Ao = 0, b) = n|, Ao / 0, c) > nj., Aq = 0, and 
d) > ni, Ao 7^ 0. The thin dashed lines indicate the 
particle dispersions below the Fermi level ep- 

Note that, when the densities in both spin states are 
equal (corresponding to Sep = 0), the dispersion rela- 
tion reduces to the usual Bogoliubov dispersion hcui^ — 
y/S,^ + |Aop describing particles above the Fermi level, 
i.e. Ek > and holes (for which the dispersion is given 
by minus the particle dispersion) below ep. It is clear 
that the Bogoliubov transformation couples particles in 
state I a) with holes in state |— a), see for example Fig. 
and that for unequal densities the dispersion relations are 
shifted with a constant ±6€p/2 such that there appear 



two separate branches in the excitation spectrum of the 
Bogoliubov quasi particles as shown in Figs. |^c and ^d. 
For n| > n^, the negative sign of TtWkj around the Fermi 
level e p indicates that the energy states are partially filled 
with spin down holes below ep, and with spin up elec- 
trons in a small region above the Fermi level. Therefore, 
the lower branch is gapless when Ag < Sep/2, whereas 
the upper one always has a gap, even at Ag = 0. The 
case of unequal densities is thus analogous to a gapless 
superconductor. 

By plugging the transformation Eq. ( p^ into Eq. 
it is easily verified that the densities satisfy 

na = ^ XI {l"kP^(?i'^k,a) + IwkHl - Nihuji,^.^)]} , 

k 

(15) 

where N{huJk^a) = l/(exp [/J^t^k^a] + 1) = (^k,a^k,a) is 
the Fermi distribution for the Bogoliubov quasi particles, 
and (3 — l/ksT. For fixed Ua, Eq. ( p5|) determines the 
chemical potentials fi'^ of the particles in state \a). 

Subsequently, the equilibrium value of the BCS order 
parameter is calculated from Aq — Vq {tpi (x)?/;| (x)) . Sub- 
stituting Eqs. ( |l0|) and (|l3| ) for ip^,i{x.), this leads to the 
following BCS 'gap equation' 



V ^ 



1 - N{riuj^^^) - N{hw^,x) 



This equation has an ultra-violet divergence, as a conse- 
quence of the fact that we made the assumption that the 
interparticle interaction is local, i.e. momentum indepen- 
dent. However, from the Lippmann-Schwinger equation 
for the two-body scattering matrix [1L9[ 



1 



1 



-y— 



V 



(17) 



we find that this divergence is cancelled by a renormaliza- 
tion of 1/Vb to 1/T^^ [|l^, and the gap equation becomes 




W 2ek 



1 



rp2B ■ 

(18) 



Eliminating from this equation both chemical potentials 
/i„ by means of Eq. (|l5|), and equating Aq to zero, one 
finds the critical temperature Tc as a function of both 
hyperfine densities in the gas. If the hyperfine densities 
are taken to be equal, the critical temperature can be 
calculated analytically |20|, resulting in 



Tr. 



Sep 



exp 



2kp\a\ 



(19) 



where 7 = 0.5772 is Euler's constant and kp ~ ^2mep/h 
is again the wave vector corresponding to the Fermi en- 
ergy ep. Including fluctuations changes only the prefac- 
tor of Eq. (|l^) g. Although this is expected to lower 
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the critical temperature somewhat, the exponential de- 
pendence of Tc on the scattering length is most important 
for our purposes. Since taking fluctuations into account 
self-consistently is rather difficult, in particular in the 
inhomogeneous case, we will come back to the effect of 
fluctuations on the transition to a superfluid state in a 
future publication and consider here only the mean-field 
theory, which is also known as the many-body T-matrix 
theory. 

As mentioned previously, the densities of particles, and 
hence the chemical potentials, need not be equal in both 
spinstates. In Fig. ^ we plot several contour plots of 
the critical temperature for the homogeneous gas in the 
n| — ni plane. As can be seen from this figure, the most 
favorable situation is that, given a certain total density 
of atoms, both hyperfine states are equally occupied be- 
cause this gives rise to the highest critical temperature. 
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FIG. 4. Contours of the critical temperature as a function 
of the hyperfine densities and n^, for X) T = 0.01 nK, 2) 
T = 11 nK, 3) r = 37 nK, and 4) T = 1725 nK. The dashed 
line is the spinodal line. 

When the two hyperfine states are not equally occu- 
pied, it can be shown that there is a nonzero critical tem- 
perature only when the spin 'polarization' |n| — nj/(n| + 
ni) < 3kBTc/2ep. Also, for fixed average Fermi level 
ep and increasing difference Sep, the critical tempera- 
ture decreases, and there is no transition at all when 

Sep ^ /cbTc(O), with Tc{0) the critical temperature when 
Sep — 1^. This behavior is similar to what occurs in 
superconductors placed in a magnetic field and can be 
understood physically from the fact that the formation 
of Cooper pairs spreads the occupation of energy levels 
only over an energy interval of order Aq ~ fc^Tc around 
the respective Fermi levels /x'^ and ii'^. Moreover, pairing 
between atoms at the average Fermi energy can only take 
place if there exists an overlap between the Fermi distri- 
butions of the two spin states in this region of momentum 
space. This indeed shows that Sep must be smaller than 
about ksTciO). 



The dashed line in Fig. ^ is the spinodal line, above 
which the gas becomes mechanically unstable. We will 
return to this issue in the next section. 

C. Mechanical stability of a two-component Fermi 

gas 

As was already pointed out in Ref. an important 
requirement for a BCS transition to occur, is that the sys- 
tem is mechanically stable against density fluctuations. 
The negative s-wave scattering length induces an effec- 
tively attractive interatomic potential, so if the density 
of particles becomes too large, the system can collapse 
to a fluid or solid state before the systems becomes su- 
perfluid in the (metastable) gaseous phase. In general, 
for mechanical stability of the gas at the critical tempera- 
ture, we must require that the velocities of the two sound 
modes in the normal state of the gas are real. These ve- 
locities can be calculated from the free energy density / 
of the gas. Since the temperatures of interest are so low 
that ksT <^ ep, we can consider the zero-temperature 
limit, in which the free energy density amounts to the 
average energy density / = {E)/V. We thus have 

f = fa + fint 

— /r 2\2/3/ 5/3 , 5/3x ^ , rr2B /on^ 

^10^^^ K )~+"'T"-ir ' (20) 

where /o is the ideal gas free energy density of the parti- 
cles in each hyperflne state at T = 0, and fint is the free 
energy density that arises due to interactions between 
particles in both spin states. The corresponding sound 
velocities squared are determined by the eigenvalues of 
the matrix 

/ 9"/ \ 

1 d\f d\f I ' 

\ dnidn-^ dnydn^ / 

leading to the condition that nin^a^ < (7r/48)^. The line 
in the n-f — Uj plane, where the equality holds, is called 
the spinodal line, and for the homogeneous ®Li gas it is 
plotted as the dashed line in Fig. ^. 

Notice, however, that a spin-polarized Fermi gas be- 
comes unstable at densities above the spinodal line, ir- 
respective of the sign of the scattering length a. There- 
fore the question arises as to what exactly happens at 
densities above the spinodal line, and whether there is a 
difference in the behavior for positive or negative s-wave 
scattering length. First of all, notice that the matrix 

f / duadnp has an eigenvalue A = at the spinodal 
point. The corresponding eigenvector eo points into the 
unstable direction of the phase space. For equal den- 
sities of the two hyperflne states, it is straightforward 
to calculate that eo = 1/a/2(=f1, 1), where the upper and 
lower sign refer to positive and negative scattering length 
a, respectively. We therefore conclude that for a nega- 
tive s-wave scattering length, the gas collapses to a dense 
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phase (probably a solid), whereas for positive a it phase 
separates into two dilute gaseous phases with opposite 
'magnetization'. Since the last situation might be of in- 
terest for future experiments with other fermionic atoms 
than ^Li, we consider now for a moment also the a > 
case. 



1. The a > Q case 

To analyze the stability at positive a, we notice that 
the pressure of the gas at zero temperature is given by 
p = -d{E)/dV. We thus find that 



P=PO+ Pint 

= \i6^)V^inf + nf)^+n,n,T'''. (21) 



Introducing for future convenience dimensionless vari- 
ables according to a; = n^a'^, y = n^a^, Af|^j = 
(2TOaV?i^)AiTa> P = a^{'2ma'^/n'^)p, and F = 
a^{2ma'^ /h^)f, it follows from Eqs. ^ and (|^ that 



Pi^X, y) = -(6^2)2/3(^^5/3 ^ ^5/3) ^ 8^^^ (22) 

5 

F{x, y) = -(6^2)2/3(^5/3 ^ y5/3^ ^ 8^^^ (23) 
5 



M^{x, y) = (67r2)2/3a;2/3 _^ g^^y 
Mx(a;,y) = (67r2)2/3y2/3^8^^^ 



(24) 
(25) 



where we used that — df /dn^^i. Notice that these 
equations are symmetric under the exchange of the vari- 
ables X and y, or rather the index "f and J,. 

The condition that must be fulfilled for a phase sep- 
aration, is that an unstable phase U separates into two 
distinct phases 5*1 and S2 in the stable region of the phase 
space, in such a way that both the pressure as well as the 
chemical potential in the two stable phases are equal. 
Since in our case we are dealing with a gas consisting of 
two constituents, we require that both chemical poten- 
tials and must be equal in the two stable phases, 
otherwise particles would still prefer one phase above the 
other, and there would be no equilibrium. A third condi- 
tion that must hold is that the total number of particles 
in each spin state must be conserved. In Fig. ^ we show 
the spinodal line, in terms of the dimensionless variables 
x,y, i.e. xy = (jr/AS)'^. Furthermore we plotted an un- 
stable point U, which separates into points Si = {xi,yi) 
and 5*2 = (2:2,2/2) in the stable regime of phase space. 
Next we will deduce the exact position of these points 
Si and 5*2 from the above mentioned conditions on the 
phase separation. 

From the condition on the pressure and the symme- 
try of Eq. (H) it follows that Pg, = P{xi,yi) = Ps^ = 
P{x2,y2) is satisfied if xi — 2/2 and X2 = yi- In other 
words, the separation points Si and S2 lie symmetric 
in the 7i| — ni plane. The condition on the chemical 



potentials, i.e. M|^|(xi,i/i) = A^T,i(-^2, 2/2) now deter- 
mines the exact position of the points 6*1 = (a;i,yi) and 
'S'2 = (2^27 2/2 )• From the symmetry of Eqs. ( p^ ) and 
(H) we see that M|(xi,yi) = Mi{yi,xi) = (0:2,^2), 
so M^{xi,yi) = Mi{xi,yi) in each individual separation 
point Si, i — 1,2. This shows that, in practice we 
are looking for intersections of the curves M|(a;,y) = 
Mi{x,y) = M. Again from symmetry, it is immediately 
clear that there is always a point of intersection of the 
two curves somewhere on the line x — y, but for certain 
values of M there can be two additional points of intersec- 
tion, which are plotted as the dashed line in Fig. |[ This 
line is the phase-separation line. As we will prove later 
on, it coincides with the spinodal line at a; = y = 7r/48, 
and lies below the spinodal line, in the stable region of 
phase space elsewhere. 



Phase-separation line 
Spinodal line 




FIG. 5. Plot of phase-separation line (dashed) as a func- 
tion of the dimensionless densities x = n^a'^ and y = rija"^, 
together with the spinodal line xy = (7r/48)^, above which 
the gas phase separates to the dashed line. As an example, 
the unstable phase U separates to the stable phases Si and 
S'2, with volume fractions Vi = riV and V2 ~ {l — 'q)V, respec- 
tively. Note that in the regions between the phase-separation 
line and the spinodal line, the gas is metastable. 

The third condition requiring conservation of the to- 
tal number of particles in each spin state determines the 
volume fractions Vi/V and V2/V of the two phases. For 
an unstable homogeneous system of volume V, and with 



V and TV 



V particles in the two hyperfine 



states, we have that after the phase separation 



= n'^^Vi+n^W2 



nf'Vi + 



1^2. 



Of course, the total density is also constant so we have 

T r TT Q Q CO 

n-totai = ""t + "'X = iT'-^^ + = n-^^ + n^^ , which means 
that the points U, Si and S'2 must lie on a straight line 
given by + ni — ntotah as is indicated for the points 
U, Si, and S'2 by the dotted hne in Fig. ||. Defining 
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^2/- ■ and 



otal 1 



now = /ntotaU = n^i^ /ntotai - "| 
r] = {P^ - I3^)/{1 - 2(3^), we find after a little algebra 
that Vi = TyF, and V2 = (1 — ?7)^^- So, the phase separa- 
tion is such that for arbitrary position of the point U on 
the unstable part of the dotted Hue in Fig. |[ the system 
separates to the same two stable points and 5*2; the 
exact position of U determines only the volume fractions 
of the stable phases. The phase points Si and S2 have 
the same total density but differ in 'spin magnetization' 
by an amount \n^'^ — n^^^\. Therefore, the phase separa- 
tion corresponds to a spin decomposition that is driven 
by the fact that at sufficiently high densities the loss in 
interaction energy between the two species compensates 
for the gain in kinetic energy due to the Pauli exclusion 
principle. 

To gain even more understanding in this phase sep- 
aration, and to distinguish later on the situation with 
negative a from the case with positive a, we consider the 
dimensionless free energy in Eq. (23) more closely. It is 
clear from Fig. ^ that the phase separation takes place 
on lines x -\- y = constant. Therefore we introduce new 
variables n and z such that 

X = n — z 
y = n + z, 

i.e. the n-axis lies along the line x = y in Fig. ^, and the 
z-axis lies along the line y — —x. Lines of constant x + y 
therefore have a constant n (density), and run parallel 
to the z-axis. The dimensionless free energy F{x, y) in 
terms of these new variables now becomes 



F(n, Z) = 1(6^2)3/2 _ ^)5/3 _^ _^ ^)5/3 

5 L 



(26) 



Note that, since the original variables x and y must 
be positive, also ?t. > 0, and for given n, we have 
—n < z < +n. Taking the derivative of F{n, z) with 
respect to z at constant n, it is found that 



From Fig. |6| we see that the maximum at z = for 
fixed n > risp is flanked by two minima in the free energy, 
which move outward in the ±z-direction for increasing n. 
Moreover, for n = ric = 97r/256 the minima just appear 
at z = ±n, i.e. at the y-axis in y = 97r/128 and at the 
X-axis in a; = 97r/128, respectively, in the original dimen- 
sionless density variables x and y. The important point is 
now that these two minima in the free energy F{n, z) for 
fixed n are, after transforming back to a: — y coordinates, 
precisely the stable separation points S'l and S2- Because 
of symmetry, they obey all conditions that we imposed 
on them. Furthermore we notice that for n > ric, or 
total density ntotai > 97r/128a3, the spin separation is 
complete, i.e. one part of the volume is occupied with 
only atoms in hyperfine level | t), the rest of the volume 
only contains atoms in state | |). The densities of both 
phases is in this case evidently n 
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FIG. 6. Plots of the dimensionless free energy F{n, z) as 
a function of z = (y — for a) n = 7r/96 < risp, 

b) n — 7r/48 — Uap, c) n^p < n = 67r/256 < ric, and 
d) n — 97r/256 = ric- 



oz 



[n - zf/'^ + {n + z) 



2/3 



— IGttz, 



which is zero at z = for all values of n. Hence, there is 
always an extremum in the free energy F(n, z) at the line 
z = 0. To see whether this is a minimum or a maximum, 
we have to analyze the second derivative 



d'^F 



dz^ 



z=0 



,1/3 



I67r, 



which is positive for n < Usp — 7r/48, zero at n — Usp, 
and negative for n > Usp. So, the minimum in the free 
energy F(n, z) at constant n and z = changes into a 
maximum at n = Ugp, which exactly coincides with the 
spinodal point at x — y. This behavior is shown in Fig. 
where we plot F{n, z) for a) n < Usp, b) n = rigp, c) and 
d) n > Usp, as a function of z. 



2. The a <0 case 

We now consider the case where the scattering length 
a < 0, as is the case for the ^Li system. Introducing 
again dimensionless variables according to x = J^tI^P 



and y — ni\a\ , and after the substitution x 



and 



y = n + z, respectively, the dimensionless free energy is 
readily seen to be 

Z) = -(6^2)3/2 [(^ _ ^)5/3 ^ ^ ^)5/3' 

5 L J 

-8TT{n^ - z^). (27) 

The first derivative of F in the z-direction is given by 
OF 



{n - z)^/^ + {n + z) 



2/3 



IGttz, 
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which is always zero at z = 0. The second derivative 
with respect to the variable z at z = is given by 



2\2/3 



,1/3 



167r, 



which is for all allowed values of n larger than zero. 
Therefore we conclude that there indeed can be no phase 
separation in the z-direction along the lines n ~ constant 
as was the case for positive a. 

Instead, the phase separation in the unstable region of 
phase space above the spinodal line takes place in the n- 
direction. This can be shown by considering the second 
derivative of F with respect to n, i.e. 



1 



1 



{n - z)i/3 + zy/3 



- 167r 



which at z = 0, or a; = y, becomes zero exactly at 
n = Usp = 7r/48. The fact that the second derivative of 
the free energy is zero at some point signals an instability 
in that direction (in the a > case, the second deriva- 
tive of F with respect to z just became zero at n = n^p). 
So we find that in the case of negative scattering length, 
the unstable point U in phase space will separate into a 
phase 5*1 with lower total particle density, and a phase 
^2 with higher total particle density, or in other words, 
to a gaseous and a dense (solid) state. However, we do 
not have an appropriate theory that can also describe 
the dense phase. Therefore we do not consider this kind 
of phase separation, which is very common in gases and 
liquids, further here. 



III. INHOMOGENEOUS FERMI GAS 

A. Local density approximation 

Until now we only considered a homogeneous gas of 
spin-polarized atomic ^Li. In reality, however, exper- 
iments with ultra-cold atomic gases are performed by 
trapping and evaporatively cooling the gas in an external 
potential which generally can be modeled by an isotropic 
harmonic oscillator V{r) = imw^r^, where w is the trap- 
ping frequency. An exact calculation of the (inhomoge- 
neous) density of the gas by calculating all eigenstates of 
the trapping potential is very elaborate but has neverthe- 
less been performed for the bosonic isotopes ^Li [ pT|j22| 
and ^^Rb [Q. Fortunately, in the fermionic system it is 
a good approximation to make use of the local density 
approximation, which treats the system as being locally 
homogeneous. This requires in the first place that the 
correlation length ^ = 0(l/fci?) is much shorter than 
the length scale / = y^h/muj over which the densities 
change. This condition is equal to the condition that 
the level spacing huj of the trapping potential is much 
smaller than the Fermi energy. Secondly, below the crit- 
ical temperature, the size of the Cooper pairs must be 



smaller than I, or else the trapping potential would in- 
fluence the wave function of the Cooper pairs. This size 
is essentially temperature-independent and of the order 
of 0{hvF/7TAo{Q)), where Ao(0) is the zero-temperature 
value of the BCS order parameter, and vp — hkp/m the 
Fermi velocity corresponding to e^. Of course, the local 
density approximation always breaks down at the edge of 
the gas cloud where the density vanishes and the effec- 
tive Fermi energy becomes zero, and also at the critical 
temperature where the correlation length ^ diverges. So 
at a nonzero temperature below Tc there are two spa- 
tial regions where the local density approximation is not 
valid, i.e. around the position where the local BCS order 
parameter vanishes and around the position where the 
local Fermi energy vanishes. However, these regions are 
so small, that we do not expect any important changes 
in the functional behavior of physical properties at the 
crossover from outside to inside these regions. As a re- 
sult we believe that it is rather accurate to apply the 
local density approximation to calculate Tc |^,^ . 

In this approximation, the densities n-^ and of the 
two hyperfine states together with the gap Aq can st ill be 
calculated by means of the equations derived in Sec. II B, 



with the understanding that now the effective chemical 
potentials, and consequently the densities and Ag, are 
spatially dependent through 



2B 



(28) 



where fia is the overall (constant) bare chemical potential 
of atoms in hyperfine state \a) . So, given T, /i| and /i| (or 
equivalently T, = J dr n|(r) and = J dr n|(r)), 
one can determine the values of ?^x(r), n.|(r), and Ao(r) 
self-consistently for every position r in space, as if the 
system were homogeneous. This procedure will be used 
in the next section to calculate the critical temperature 
of the spin-polarized g cls cis Oil function of the number of 
particles in the trap. 



B. Critical temperature 

The critical temperature Tc of the gas is such that 
at the center of the magnetic trap, where the density 
of the gas is highest, the energy gap Ao(0) just becomes 
nonzero for a given number of particles Nf and . First 



N. 



In Fig. 



we will consider the case where = 
the solid line shows the result of our calculation. The 
dashed line in this figure gives the critical temperature 
for the Fermi gas if one does not include the effects of the 
mean- field interaction in Eq. (|2^). In this approximation, 
the number of particles in each hyperfine state is, with 
a high degree of accuracy, given by the zero-temperature 
result iV|,| — {fi^^i/huj)^ /6, and the density in the cen- 
ter of the trap is ?^t,x(0) = (2m/i|,x/?i^)^/2/(67r2), which 
is considerably smaller than in case that the mean-field 
interaction is taken into account. As a result, the criti- 
cal temperature obtained in this manner is substantially 
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lower for an equal number of particles. From an experi- 
mental point of view, it is therefore important to include 
interactions to obtain a reliable estimate for the criti- 
cal temperature as a function of the number of trapped 
particles. 

We found that, as is also the case for a Bose gas in 
a harmonic trap [ p^ , the critical temperature, or rather 
the dimensionless parameter q/Xt^ is a universal func- 
tion of N'^/^a/l, with the thermal DeBroglie wavelength 

At = Jp.-K'h^ /mkBT and / — ^JhjmLj. The solid hue 
in Fig. [7| is can be fitted numerically very well with the 
expression 



-^—\ = 0.037exp(-1.214/-iV-s -H2.990^iV^ 
\tJ I \a\ I 

for the whole range of parameters shown in Fig. 0. 



n„(0)A|,^ = Fa 



So the central density of each spin state times the ther- 
mal wavelength is also a function of the dimensionless 
parameters Sep/ksTc and ep/ksTc- Combining these 
two equations, it follows that a/Xx^ is directly related to 
the densities in the center of the trap, i.e. 



A 



F[nt(0)A3,^,ni(0)A|,J. 



(29) 



To prove now that a/Xx^ is a function of N^^^a/l, it 
should be noticed that in general in the local density ap- 
proximation for T > Tc 

"aWAf, = /3/2 [exp{/3[/ia - n_„(r)r2S - V{r)]}] , 

(30) 
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where /3/2[2(r)] is the Fermi function originating from in- 
tegration over momenta and analogous to the Bose func- 
tion 33/2(2). Applying this equation at r = and T = 
we find that both chemical potentials are functions of 
q/Xt^ and the central densities of both hyperfine states, 
and obey 



ksTc 



= F' 



— ,nT(0)A3 ,n^(0)A3 



For a general value of r, but still at T = Tc, we can apply 
the substitution 



2kBT, 



(31) 



FIG. 7. Critical temperature as a function of the number 
of particles (solid line), when there are A'' particles present in 
both spin states. The dashed line represents the critical tem- 
perature for a gas whose density distribution is not altered by 
mean-field interactions. 

The fact that the critical temperature is a universal 
function of the parameter N^/^a/l can be understood 
easily by rewriting the gap equation Eq. ( p^ ) at the crit- 
ical temperature in the form 

\/7r Atc _ 



dxyfx 



Nj-Sep + || - ep]) + NjSep + If - ep]) 
2\x-(3ep\ 



where N{x) — l/(exp [f3x] -I- 1) is the Fermi distribution. 
This shows that at the critical temperature, a/Xx^ is a 
function of Sep/ksTc and ep/ksTc only. Equivalently, 
from the density for each spinstate given in Eq. (^5|) , and 
the fact that at the critical temperature the densities in 
the center of the trap UaiG) are critical, we find that 



in Eq. (pOj), from which it follows immediately that for 
each hyperfine state 



na(r)Af,^ 



— ,nT(0)Af,^,n|(0)Af,^,y2 



To find the total number of particles in each hyperfine 
level, we then integrate this result over the spatial extent 
of the gas cloud, resulting in 



47r / dr r na{r) 
Jo 

47r /2kBT^.^^^^ 
X^ \ mcj2 



dy { X 



F" 



a 
Atc 



n^{0)Xl^,n^{0)Xl^,y'' 



I 

XPa 



F" 



— ,n^{0)X^r^,n^{0)X^T^ 



(32) 



Multiplying Eq. ( |3^ ) on both sides with (a/l)^, and us- 
ing the result of Eq. (^9|) , it is proven that at the critical 
temperature 
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F 



1/6 a 1/6" 



(33) 



C. Superfluid state 



so that, when /i^ = /i| the dimensionless parameter a/X^^ 
is a universal function of N^^^a/l. 

The spinodal point in this case is given by N^^^a/l ~ 
0.66, and is not included in Fig. |^, because for ^Li trapped 
in a harmonic potential with frequency v — w/27r = 
144Hz, OT huj/ks — 6-9 nK, corresponding to the present 
experimental conditions at the Rice experiment Q , spin- 
odal decomposition only occurs with as many as 5.8 x 10^ 
particles. 

For an unequal number of particles in each hyperfine 
state, we find a universal surface for a/Xx^ as a function 

of N^^^a/l and N^^^a/l, as Eq. ( |3^ ) shows. However, 
since we are in this paper mainly interested in trapping 
^Li atoms, we will calculate several contours of the crit- 
ical temperature of such a gas trapped in an isotropic 
harmonic oscillator with v = 144Hz. The results are 
plotted in Fig ^. Again we see that given the total num- 
ber of particles in the gas, the most favorable situation is 
the one with equal numbers of particles in each hyperfine 
state. 




FIG. 8. Critical temperature as a function of the number of 
®Li atoms in each hyperfine state. Curve 1 to 3 give the combi- 
nations {Ni,Ni) corresponding to Tc = 3 nK, 2) Tc = 11 nK, 
and 3) Tc = 37 nK. For equal number of particles in each hy- 
perfine state, the density of particles in the center of the trap 
corresponds to 1) n|,|(0) = 0.5x10^^ cm"^, 2) 1x10^^ cm"^ 
and 3) 2 X 10^^ cm~'^ respectively. 

An important experimental question is how we could 
observe whether or not the gas is superfluid at a certain 
temperature. An immediate possibility that, in view of 
the results with the BEC experiments, comes to mind is 
to consider whether there is a change in the density pro- 
file at the critical temperature. In the next section we 
will therefore concentrate on the superfluid state of the 
gas, and determine the density profiles and in addition 
the spatial dependence of the energy gap Ao(r). 



In Fig. H the density profile n|(r) — n|(r) and the 
energy gap Ao(r) is plotted for several temperatures 
below and at the critical temperature for a gas with 
N-^ — Ni — 2.865 X 10^ particles in both hyperfine states. 
The dotted line in Fig ^ shows the density distribution 
for a gas with the same number of particles, but with 
a — instead of a = — 2160ao. It is clearly visible that 
the effect of the interaction on the density is rather large. 
Indeed, due to the attractive interactions the particles 
are pulled to the center of the trap, and the density is 
there considerably increased which is good from an ex- 
perimental point of view because it significantly increases 
the critical temperature. 
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FIG. 9. Density distribution n|(r) — n|(r) and energy gap 
Ao(r) for a ^Li atomic gas consisting of 2.865 x 10^ atoms 
in each spin state at a) T = 15nK, fo) at T = 33nK, slightly 
below Tc, and c) at T = Tc — 37nK. The left scale of each 
plot refers to the density, the right scale to energy gap. The 
open circles in b) represent Eq. (M) , and the dotted line in c) 
shows the density distribution for a gas with the same number 
of particles and at the same temperature, but with a — 0. 



Fig. hd is a similar plot, but now with unequal num- 



ber of particles in each spin state: A'^t 



3.08 X 10^ 



and Ni = 2.65 x 10 , so the total number of particles is 



the same as in the previous case. From Fig. [LOp it can 
be seen that the presence of the order parameter tends 
to decrease the difference in densities of each hyperfine 
level. This can physically be understood from the fact 
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that the most favorable condition for the formation of 
Cooper pairs is that both densities are equal. 
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FIG. 10. Density distribution n^{r), n|(r) and energy gap 
Ao(r) for a ®Li atomic gas with N-f = 3.08 x 10^ and 
Ni — 2.65 X lO'' below and at the critical temperature. In 
a) r = 15 nK, fe) at T = 30 uK, slightly below Tc, and 
c) at T = Tc = 34 nK. The left scale of each plot refers 
to the density, the right scale to energy gap. The dotted 
lines in c) represent the density profiles for a non-interacting 
gas, with the same number of particles A''| = 3.08 x 10^, and 
Ni = 2.65 X lO'^, respectively. 

The most important observation that we can make 
from both Figs. || and is that there is almost no change 
in the density of the gas going from the normal to the 
superfluid phase. This also leads, as will be explained in 
more detail below, to the conclusion that a measurement 
of collective excitations will not give a good signature for 
the presence of a superfluid state |26j . A second obser- 
vation is that from the BCS theory in superconductors 
| p6[ , it is known that the order parameter Aq close to the 
critical temperature vanishes as 



Ao(r) ~ 1.74Ao(0)Vl-T/re, 



(34) 



where Ao(0) is the zero-temperature value of Aq, which 
in turn is related to the critical temperature as 



Ao(0) ~ 1.76fcsr, 



(35) 



For Fig. |9|a it follows from Eq. (19) that the critical 
temperature corresponding to the density of the gas in 
the center of the trap is much larger than the temper- 
ature (T = 15 nK) itself. Hence the value of the or- 
der parameter approaches the zero-temperature limit in 



this case. Using that rc['^(0)] ~ 37nK, one finds from 
Eq. (H) that Ao(0)/fcs = 65.1nK. Comparing this with 
Ao(r = 0)/kB = 65.0nK, we find that there is indeed 
a rather good agreement. Also, since the temperature 
in Fig. ^ is only slightly below the critical tempera- 
ture, we can compare relation Eq. (^^ with the val- 
ues of Ao(r) in this figure. At each spatial position 
where Ao(r) > 0, we can from the local value of the 
Fermi energy ei?(r), extract the local critical tempera- 
ture Tc{r) from Eq. ( |l9|) and use Eq. ( |35| ) to compare 
Ao{T\r])/kB = 3.06rc(r)^l-T/rc(r) (open circles in 
Fig. !&) with Ao(r). Again, the agreement is very good. 

Finally, we want to check that our local density ap- 
proximation is indeed valid under the conditions of in- 
terest. From Fig. ^, one finds that at r = 0, the value of 
1/kpl ~ 0.06, and for example at r = 0.05 mm, we find 
l/kpl ~ 0.64. Furthermore, from the zero-temperature 
value of Aq in Eq. (]35|), we find that the size of the Cooper 
pairs relative to the trapping parameter I is about 0.58, 
so the local density approximation starts to break down 
if we are far below the critical temperature. In that case 
a more accurate approach is required, at least for the 
relatively large trapping frequencies used here. 



IV. DISCUSSION AND CONCLUSION 

As mentioned above, an important experimental prob- 
lem is the detection of the superfiuid state. In contrast 
to the Bose- Einstein condensation experiments, there is 
no clear signature in the density distribu tion w hen the 
gas becomes superfiuid, as shown in Sec. II 
fore, a measurement of the collective excitations, 



There- 
or den- 
sity fluctuations, will not provide useful information on 
the presence of the superfluid phase as well. This can 
also be understood from the dissipationless (linear) hy- 
drodynamic equations governing the density fluctuations 
in the system. Considering only the optimal situation of 
an equal number of particles N-^- = N^, these equations 
are given, for a gas trapped in an external potential V{r), 
by 



d\" 1 

+ — (Vp-f nVy(r)) = 
at m 

ot m 



(36a) 
(36b) 
(36c) 
(36d) 



where n = n„ -|- ris is the total density that consists of 
a normal and superfluid part, j" = nnV„ -I- n^v^ is the 
density current with (v„) the superfluid (normal) ve- 
locity, e the average energy density, = fip + TsVn the 
energy current, and s is the entropy density p^,E5|. Note 
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that the same equations in fact also describe the collec- 
tive modes of a trapped Bose-condensed gas |28 2^. 

To show that these hydrodynamic equations result in 
identical equations for the collective excitations in the 
normal and the superfluid phase, we first of all note that 
for the densities of interest the gas will be in the hy- 
drodynamic limit, meaning that the time scales for the 
density fluctuations (which are of the order of the in- 
verse trapping frequency) are much slower than the time 
between elastic collisions. In this hydrodynamc regime, 
density fluctuations and temperature fluctuations influ- 
ence each other with a coupling proportional to cy — Cp, 
where cy (cp) is the heat capacity per particle at con- 
stant volume (pressure). However, for the very low tem- 
peratures of interest, one can assume that the heat ca- 
pacities of the gas satisfy cy ~ Cp. Indeed, for a ho- 
mogeneous Fermi gas, one finds in the limit T I that 
(cp — cv)/cv = OdkBT/ep]'^) and is thus very small. 
As a result, the density and temperature fluctuations are 
effectively uncoupled [2^. As a consequence, Eq. (36c), 
descibing second sound, decouples from the other 3 equa- 
tions, and it suffices to consider only density fluctuations 
n{r,t) — no{r) + Sn{r,t). Note also that if we have an 
unequal number of particles, i.e. 7V| 7^ N^, the density 
fluctuations are coupled to fluctuations in the 'magneti- 
zation' n| — nj^ and we need to generalize these equations. 
For equal number of particles thes e 'sp in waves' decouple, 
however, as we have seen in Sec. [IC. 



In the normal phase, the Josephs on r elation Eq. (36d) 
must be dropped. Linearizing Eq. (36t), which is in fact 
just Newton's law, we arrive at 

, , 9v(r, t) 1 „ / r , , , dp ^ , , 
at m \ on 

-— no(r)VT/(r) - — (5n(r, t)Vy(r) 
m m 

m \on / m 

where in the second line we used that in equilibrium 

Vp[no(r)] = -no(r)Vy(r) (37) 

and furthermore that the pressure p is a function of the 
density only at fixed temperature. This result, together 
with the continuity equation in Eq. ( |36a| ) describes first 
sound in a trapped Fermi gas. 

In the superfluid phase, first sound has Vs = v„ and 
Eq. (36b) now becomes 



no{r) 



dt 



-vfen(M) 
m \on 

m V on 



-Sn{r, i)Vy(r) 



VV{r) 6n{r,t) 



--^V6nir,t). 
m on 



(38) 



However, in the s uperfl uid phase we also have the Joseph- 
son relation Eq. (36d) for the superfluid velocity. Using 



that the local chemical potential /i = fio[nQ{r)+Sn{r, t)] + 
V{r), where /xg is the homogeneous chemical potential 
including the effects of interactions, and linearizing also 
this equation leads to 



"o(i")^ = V ^lo[no{v)] + —Sn{r,t) 

at m \ on 



^ 5n(r,i)V— — + — — V(5n(r,i) 

m V on on 



y(r; 



(39) 



because in equilibrium the chemical potential must sat- 
isfy 



/io["o(r)] + ^(r) — constant. 



Since in general 



d-p 
dn 



no' 



djM)_ 
dn ' 



(40) 



(41) 



the second term on the right-hand side of Eq. (p9|) equals 
the last term on the right-hand side of Eq. (|38|). More- 
over, the first term on the right-hand side of Eq. ( ^9|) can 
be rewritten with Eq. (41) as 



^ djiQ _ dp 1 
dn 



dn nQ(r) 
1 



Vno(r) 



no(r) 



1 ^dp 

] V — 

no(r) dn 

dn 



no{r) 



V|+Vnr) 



where we used again the equilibrium condition in 
Eq. (^. We thus find that below the critical temper- 
ature, the Josephson relation Eq. (36d) is identical to 
the momentum equation (36b). As in the normal phase, 
first sound can below Tc thus be described merely by 
Eqs. (36a) and (36b). In combination with the results 



that the density profiles in the normal and superfluid 
phase are almost equal, we conclude that the hydrody- 
namic equations that describe the density fluctuations 
are almost identical, and therefore, that there will be no 
signiflcant difference in the collective excitation spectrum 
in the superfluid and normal phase respectively. Con- 
sequently, other means of experimental detection must 
be investigated. Of course, this conclusion is based on 
experiments that couple directly to density fluctuations 
such as in Refs. |^,^. If one can couple also to second 
sound one would of course observe an additional mode 
below Tc- 

Another possible way to detect the superfluid state is 
by a measurement of the two-body decay rate of the gas. 
Note that, in our case, three-body recombination pro- 
cesses are strongly suppressed, since we are dealing with 
fermions, and only have two different hyperfine states 
occupied. Above the critical temperature, the two-body 
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rate constants are essentially independent of T, and the 
magnitude is depicted as a function of the applied mag- 
netic field in Fig. ||. In analogy to the case of a Bose gas, 
where the presence of a condensate decreases the decay 
rate due to two-body processes by a factor of about 2 
p2|-^ , we now analyze the change in the decay rate due 
to the presence of Cooper pairs in the Fermi gas below 
the critical temperature. 

Using the correlator method from Ref. jS^] , it is found 
that the decay rate constant due to two-body processes 
is given by 



G{T) = G{T^)K^'^\t), 



(42) 



where the correlator 

K(2) (t) = J_ (^t (^)^t (^)^ (^)^ (^)) (43) 

equals 1 above Tc, but increases due to the nonzero expec- 
tation value {ip-^ (x)-;/'! (x)) below the critical temperature. 
Indeed, using the transformations given by Eqs. (^) and 
(p^), it is found that 

I A 12 

(2)m - 



K('){T) = 1 



T r1 ' 



where Aq/Vq is again given by the ultra-violet diverging 
expression in Eq. (16). The question now arises what 
we should use for Vq in this expression. Clearly, the 
denominator Vq should not be considered as the zero- 
momentum component of the triplet potential. Physi- 
cally, this can be understood by the fact that this value 
of Vb does not characterize the exact non-local two-body 
triplet potential in any way. It does not even reproduce 
the correct long distance behavior for the scattered wave 
function. Therefore, a first guess would be to replace Vb 
by T^^, but since the s-wave scattering length in this 
case is much larger than the effective range ry of the in- 
teraction Vt(x), it is likely that the replacement of Vq 
by T^^ underestimates the effect of the presence of the 
superfluid phase on the decay rates considerably. 

Instead, from our procedure to remove the ultra-violet 
divergence in the gap equation, we see that Vb should 
be chosen such that it satisfies the Lippmann-Schwinger 
equation 



1 



1 

Vo 



1 



1 



k<A ^ 



where a cut-off A = 0(l/ry) is introduced. Solving 
this equation for Vb, and using that the Fermi energy 
ef < Ti^K^ /2m, we find that 



Vb = r 



2B 



1 



1 - 



2aA 



(44) 



Applying this result we can now make an estimate of 
the effect of the presence of the superfluid phase. To sim- 
plify the calculation, we will calculate the effect at T = 



for a homogeneous gas with n| = = n. We then have 
that n = kp/{6TT'^). Moreover, from the BCS theory we 
know that the zero-temperature value of the order pa- 
rameter Ao(0) is related to the critical temperature by 
Ao(0) = 1.76kBTc. So 



1 



1 



1 - 



2aA 



l.TGfcsTc 37rn 
ep Akpa 



From the functional behavior of the triplet potential 
Vt (x) , we deduce that the range of the two-body interac- 
tion rv — lOOflQ. Therefore, substituting A ~ (IOOoq)"^ 
and a = — 2160ao, we find that in the case of ^Li 
atoms, Vb — 0.07T^^. Note that, as expected, the 
Fermi wavenumber hp is much smaller than the cut-off 
value A: For a density = — 10^^ cm~^, we have 
ep ~6x 10~''A;_B, resulting in kp ~ {5x lO^oo)"^, which 
is indeed much smaller than the cut-off A. Using that 



10 



12 



cm 



(see Fig. ^) it turns out 



Tc = 11 nK forn , „ ^. 

that 1 - 2aA/7r ~ 15 and the correlator X(2)(o) ~ 7, 
which is much larger than its value above the critical 
temperature. Had we used T^^ instead of Vb, the change 
in the correlator would have been only of the order of 3%. 
Even though we do not expect that the corrections to the 
decay rates are as large as the above crude argument sug- 
gests, we do believe that it might be of the order unity, 
and may be measurable. Of course, the cut-off depen- 
dence of Vb should in some way drop out of the theory 
eventually, but for this a better theory is needed, which 
takes into account the precise details of the triplet poten- 
tial and does not make use of a pseudo-potential to re- 
place it. Work along these lines is in progress, because of 
the experimental importance to have a reliable estimate 
of the changes in the relaxation rate constants. Further- 
more, note that the correlator K^'^\T) from Eq.(|3|) also 
appears in the expression of the average energy of the 
system [Q. A measurement of this quantity has been 
done for the case of Bose gases and we believe 

that also in the case of fermions, a change in the average 
energy at the critical temperature can signal the presence 
of the superfluid phase. 

In summary, we considered a gas of atomic ®Li occu- 
pying two hyperfine states trapped in a magnetic field. 
Atoms in different hyperfine levels can interact via s-wave 
scattering. Using the most up-to-date triplet potential, 
we showed that the lifetime of such a gas with a density 
of 10^^ cm~^ per hyperfine level is of the order of 1 sec- 
ond, when a magnetic bias field of 5 T is applied. At this 
density the gas becomes superfluid at a temperature of 
about 11 nK. We also investigated the mechanical sta- 
bility of a two-component Fermi gas, and showed that, 
if the two-particle interaction is repulsive, the gas is un- 
stable for spin-density fluctuations, whereas in case the 
interatomic interaction is attractive, the gas is unstable 
against density fluctuations above the spinodal line. 

Furthermore we considered the superfluid state of 
atomic ^Li, trapped in an isotropic harmonic oscillator 
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potential, in the local density approximation. We showed 
that the critical temperature a/Xx^ is a universal func- 
tion of the quantities N^^^a/l and N^^^a/l. Below the 
critical temperature there is almost no change in the den- 
sity profile of the gas cloud. Therefore we suggest that 
the presence of the superfluid state might be signaled by 
a change in the decay rates or a change in the average 
energy at the critical temperature. 
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